Award  Number:  DAMD17-01-1-0402 


TITLE:  Beam  Delivery  Verification  for  Modulated  Electron 

Radiation  Therapy  Treatment  of  Breast  Cancer 


PRINCIPAL  INVESTIGATOR:  Maxime  Ponineau 

Arthur  L.  Boyer,  Ph.D. 


CONTRACTING  ORGANIZATION:  Leland  Stanford  Jr.  University 

Stanford,  CA  94305-5401 


REPORT  DATE:  August  2003 


TYPE  OF  REPORT :  Annual  Summary 


PREPARED  FOR:  U.S.  Army  Medical  Research  and  Materiel  Command 
Fort  Detrick,  Maryland  21702-5012 


DISTRIBUTION  STATEMENT:  Approved  for  Public  Release; 

Distribution  Unlimited 


The  views,  opinions  and/or  findings  contained  in  this  report  are 
those  of  the  author (s)  and  should  not  be  construed  as  an  official 
Department  of  the  Army  position,  policy  or  decision  unless  so 
designated  by  other  documentation. 


20031121  066 


REPORT  DOCUMENTATION  PAGE 


Form  Approved 
OMB  No.  074-0188 


Public  reporting  burden  for  this  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and  maintaining 
Jhe  data  needed,  and  completing  and  reviewing  this  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information,  including  suggestions  for 
reducing  this  burden  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington,  VA  22202-4302,  and  to  the  Office  of 
Management  and  Budget,  Paperwork  Reduction  Project  (0704-01 88),  Washington,  DC  20503 


1 .  AGENCY  USE  ONLY  2.  REPORT  DATE  3.  REPORT  TYPE  AND  DATES  COVERED 

(Leave  blank)  August  2003  Annual  Summary  (1  Aug  2002  -  31  Jul  2003) 


4.  TITLE  AND  SUBTITLE  5.  FUNDING  NUMBERS 

Beam  Delivery  Verification  for  Modulated  Electron  Radiation  DAMD17-01-1-0402 
Therapy  Treatment  of  Breast  Cancer 


6.  AUTHOR(S) 

Maxime  Ponineau 
Arthur  L.  Boyer,  Ph.D. 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

Leland  Stanford  Jr.  University 
Stanford,  CA  94305-5401 


8.  PERFORMING  ORGANIZATION 
REPORT  NUMBER 


E-Mail:  Maxime  .  popineau@s tanf ord .  edu 


9.  SPONSORING  /  MONITORING 

AGENCY  NAME(S)  AND  ADDRESS(ES) 

U.S.  Army  Medical  Research  and  Materiel  Command 
Fort  Detrick,  Maryland  21702-5012 


10.  SPONSORING  /  MONITORING 
AGENCY  REPORT  NUMBER 


12a.  DISTRIBUTION  /  AVAILABILITY  STATEMENT 

Approved  for  Public  Release;  Distribution  Unlimited 


12b.  DISTRIBUTION  CODE 


13.  ABSTRACT  (Maximum  200  Words) 

Detecting  and  monitoring  cellular  and  molecular  changes  associated  with  cancer  are 
essential  to  the  verification  of  cancer  treatment.  In  this  investigation  we  show  that 
modulated  electron  radiation  therapy  can  in  principle  conform  a  region  of  high  dose  to  a 
volume  of  malignant  breast  tissue.  Currently  available  means  to  both  measure  the  malignant 
volume  and  record  progression  of  its  response  to  radiation  are  limited.  We  therefore  chose 
to  explore  another  means  to  determine  the  breast  treatment  target;  optical  tomography,  j 
Using  the  fact  that  mammalian  tissues  transmit  light  at  a  low  level  and  emit  virtually  no 
light  at  all,  optical  signatures  conferred  on  tumor  cells  by  expression  of  reporter  genes 
can  be  detected  externally  by  photon-detecting  systems. 

In  this  paper  we  study  a  first  approach  to  a  gradient-based  iterative  procedure  using  a 
finite-difference  scheme  as  forward  model.  We  provide  results  at  the  current  state  of 
research  and  discuss  ways  to  address  current  limitations. 


14.  SUBJECT  TERMS 

Modulated  Electron  Radiation  Therapy 


17.  SECURITY  CLASSIFICATION 
OF  REPORT 

Unclassified 


NSN  7540-01-280-5500 


18.  SECURITY  CLASSIFICATION 
OF  THIS  PAGE 

Unclassified 


19.  SECURITY  CLASSIFICATION 
OF  ABSTRACT 

Unclassified 


15.  NUMBER  OF  PAGES 

10 

16.  PRICE  CODE 


20.  LIMITATION  OF  ABSTRACT 


Unlimited 


Standard  Form  298  (Rev.  2-89) 

Prescribed  by  ANSI  Std.  239-18 
298-102 


Table  of  Contents 


Cover . ....1 

SF  298 . . . . . 2 

Table  of  Contents . . . 3 

Introduction . . . 4 

Body . . . 4 

Key  Research  Accomplishments . 8 

Reportable  Outcomes . . 

Conclusions . . . 9 

References.. . 10 


Appendices 


none 


Introduction 


Detecting  and  monitoring  cellular  and  molecular  changes  associated  with  cancer  are  essential  to 
our  understanding  of  neoplasia  and  to  the  verification  of  its  treatment.  In  this  investigation  we 
show  that  modulated  electron  radiation  therapy  can  in  principle  conform  a  region  of  high  dose  to 
a  volume  of  malignant  breast  tissue.  Currently  available  means  to  both  measure  the  malignant 
volume  and  record  progression  of  its  response  to  radiation  are  limited.  We  therefore  chose  to 
explore  another  means  to  determine  the  breast  treatment  target;  optical  tomography.  Using  the 
fact  that  mammalian  tissues  transmit  light  at  a  low  level  and  emit  virtually  no  light  at  all,  optical 
signatures  conferred  on  tumor  cells  by  expression  of  reporter  genes  can  be  detected  externally  by 
photon-detecting  systems. 

The  realization  of  such  in  vivo  optical  imaging  therefore  requires  three  components:  (1)  a  fast 
and  accurate  measurement  means  of  near-infrared  (NIR)  light;  (2)  a  computationally-efficient 
reconstruction  algorithm  and  (3)  the  existence  of  a  significant  contrast  in  the  optical  properties  of 
the  medium  that  we  are  reconstructing.  For  (3),  expression  of  reporter  genes  such  as  the 
bioluminescent  enzyme  firefly  luciferase  can  effectively  modify  tumor  cells  optical  properties, 
thus  securing  sufficient  contrast  between  the  tumor  and  the  rest  of  the  medium. 

As  far  as  the  inversion  algorithm  (2)  is  concerned,  most  current  efforts  have  been  hampered  by 
the  computationally-intensive  nature  of  a  full  three-dimensional  data  inversion. 

In  order  to  achieve  a  first  simple  and  robust  method,  a  gradient-based  reconstruction  algorithm 
was  investigated  to  invert  a  bioluminescence  problem  and  compute  the  shape  of  a 
bioluminescent  tumor.  A  similar  approach  is  described  in  article  [1]. 

This  method  was  implemented  through  an  iterative  procedure  that  falls  into  four  distinctive 
steps:  we  successively  (1)  chose  an  initial  source  distribution;  (2)  computed  its  response  a 
forward  model;  (3)  computed  an  objective  function  (which  estimates  the  distance  between  our 
‘guessed’  source  distribution  and  the  actual  one)  and  computed  the  gradient  of  the  objective 
function  with  respect  to  the  source  intensity  distribution  and  (4)  moved  along  the  gradient 
direction  until  we  found  a  closer  distribution  and  updated  the  source  distribution.  The  iterative 
feature  naturally  came  from  the  fact  that  steps  (2)-(4)  was  performed  repeatedly  until  our 
computed  source  distribution  produced  a  signal  that  was  ‘close  enough’  (we  will  discuss  later 
what  that  means)  to  the  actual  measurement.  This  procedure  was  implemented  on  Matlab  and  is 
represented  in  figure  1. 
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fig  1:  flow  diagram  of  the  gradient-based  iterative  procedure 


In  highly  scattering  media  such  as  biological  tissue,  light  propagation  is  well  approximated  by 
the  photon  diffusion  equation,  which  is  as  follows: 

V2U(r,t)-^U(r,t)  =~S{r,t) ,  (1) 

D  D  dt  D 

Where  U  (r,t)  is  the  photon  fluence  rate  [photons/cm2s],  jua  and  jus  are  the  absorption  and 
scattering  coefficients,  respectively.  D  =  v/3 jus  is  the  photon  diffusion  coefficient,  and  S(r,t) is 

the  number  of  photon  emitted  by  the  source  (tumor)  at  position  r  and  time  t. 

Equation  (1)  is  actually  a  particular  case  of  a  more  general  one: 

dU_dfdU)  d(^du)  TT  c 

dt  dx  J  +  3y[D  dy )  CMa  +  ’  (  } 

in  which  we  now  assume  that  the  diffusion  coefficient  D  is  constant  and  uniform. 

At  this  point,  we  need  to  use  a  forward  model  to  compute  the  system's  response  to  a  given  source 
distribution. 


I  -  The  forward  model  and  Predicted  Measurement: 


1.  Using  the  method  described  in  article  [2],  the  propagation  of  light  intensity  is  assumed  to 
be  a  simple  spherical  damped  wave,  i.e.  described  by: 

d,(r  rV-^o  exp(4~~rj) 

D  4a\r-rs\ 

With  such  assumption,  we  can  compute  the  system’s  response  to  any  distribution  of 
source  by  just  considering  the  distribution  as  the  sum  of  local  source  and  using  the 
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superposition  principle.  Indeed,  since  equation  (1)  is  linear  in  the  source  and^(r,t),  the 
response  to  a  sum  of  point  sources  is  simply  the  sum  of  the  responses  to  each  point 
source. 

2.  The  second  approach  is  much  more  realistic  and  uses  a  finite-difference  scheme  to 
compute  the  response.  Thus  equation  (1)  is  solved  by  replacing  the  temporal  and  spatial 
derivatives  by  their  finite-difference  approximations: 

du  -ir-u*-1 

dt  At 

We  similarly  have  that: 

d(cdu)~cu‘«-<-u*>  alr3e/l  ptAn zEiti 

dx )  Ax2  '  3yv  J  Ay2 

As  explained  in  article  [1],  different  possibilities  arise  when  it  comes  to  discretizing  a 
continuous  differential  equation.  In  particular,  two  finite-difference  schemes  are 
particularly  natural:  the  explicit  scheme,  which  expresses  the  current  state  (n+1)  in  terms 
of  the  previous  state  (n),  and  the  implicit  scheme,  which  expresses  the  future  state  in 
terms  of  the  current  one.  The  explicit  scheme  is  the  simplest  since  it  computes  the  current 
state  only  from  parameters  that  we  already  know  (i.e.  the  past  state),  but  further  study 
proves  that  it  is  only  conditionally-stable.  On  the  other  hand,  the  implicit  scheme 
involves  more  computation  but  turns  out  to  be  unconditionally-stable.  I  followed  the 
method  chosen  in  Article  [1],  which  tries  to  get  the  best  of  both  worlds  by  using  an 
Alternation-Direction  Implicit  (ADI)  scheme.  In  the  ADI  method,  the  computation  of 
Un+1  from  U"  is  split  in  two  steps.  In  the  first  half-time  step,  the  spatial  derivative  in 
only  one  (say  x)  direction  is  evaluated  at  the  present  time  step  (implicit)  and  the  other 
spatial  derivative  (say  y)  is  evaluated  in  the  previous  time-step  (explicit).  In  the  next  half¬ 
time  step,  the  implicit  and  explicit  directions  are  switched. 

According  to  this  scheme,  equation  (1)  is  discretized  and  it  can  then  be  shown  that 

equation  the  relationship  between  Un  and  un+V2  can  be  expressed  in  the  following  matrix 
form: 

AJjn+ll2  =Bun  +Sn+V  2 

Un+V2  =A~1BUn +A~lSn+V2  ^ 

where  Un+V2is  an  intermediary  state.  A  similar  relationship  is  then  found  between 
Un+m and  Un+1. 

II  -  Analysis  Scheme: 

The  analysis  scheme  only  consists  in  the  calculation  of  the  objective  function  <f> ,  which  is  simply 

the  sum  of  the  squared  differences  between  the  computed  signal  and  the  actual,  measured  one. 

This  is  naturally  a  function  of  the  distribution  of  sources,  denoted  C  • 

«K<r)=I iv-u'M  <4) 

seM 
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Ill  -  Updating  Scheme: 


This  is  where  we  turn  to  the  computation  of  the  gradient.  This  part  is  naturally  the  most 
computationally-heavy.  The  derivative  of  the  objective  function  with  respect  to  the  source 
distribution  £  is  given  by: 


d(j) 

dC 


y  y  du„ 

dCr 


where  n  represents  the  time  steps  and  p  is  the  index  of  a  grid  point.  If  we  consider  a  single  time 
step,  this  simplifies  to: 


d(j) 


d</> 


dUp 

n 


(5) 


9Cr  ptW  BCr  ’ 

The  first  part,  which  is  the  derivative  of  the  error  function  with  respect  to  the  intensity,  comes 
easily  if  we  remember  that  the  error  function  only  depends  on  the  intensity  at  grid  point  p  if  p  is 
on  the  boundary,  ie  if  it  is  part  of  the  measurement  set.  In  this  particular  case,  the  error  function 
is  a  simple  quadratic  function  of  the  intensity,  and  we  thus  have: 


d* ..  <>i  \-jyr-vr) 


dUp  dUp 


0 


,pe  M 
otherwise 


(6) 


Now,  the  signal’s  derivative  with  respect  to  the  source  distribution  naturally  depends  on  the 
assumption  made  for  the  forward  model. 


1.  In  the  framework  of  article  [2],  this  derivative  follows  directly  from  our  assumption: 


d</> 

K, 


=  '£-2{yp-u') 


v  expl 


(4>~r'P 


peM 


D 


4 n\rp  -  rr 


(11) 


2. 


In  the  finite-difference  forward  model,  this  gradient  computation  is  naturally  more 
involved  and  requires  additional  care.  In  particular,  one  interesting  way  is  to  use  equation 
(3)  to  see  that,  if  we  consider  (3)  a  time  n-1/2  and  n  instead  of  n  and  n+1/2,  we  get: 


thus 


dA 

dC 

dun 

ZCr 


+  A 


=  A~ 


dUn  dB 


dCr  BCr 


u 


n-1/2 


dB 


U 


n-1/2 


?A  ^ 

dA-un 


(7) 


vdf,'  K  / 

which  can  be  easily  calculated  by  since  we  know  the  intensity  U  and  the  matrices  A  and 
B.  This  method  is  yet  to  be  implemented  this  method  but  it  promises  to  be  more  adapted 
to  our  problem  than  1,  which  is  over-simplifying. 


IV  -  Line  minimization: 

Once  we  have  the  vector  (d<p/d£r  )rgiW  for  a  current  position  £0 ,  we  can  move  along  the  direction 

defined  by  (d<f>fd£r)reM  and  look  for  the  position  that  yields  the  smallest  objective  function 

possible.  This  line  minimization  that  we  used  does  not  systematically  yield  the  global  minimum 
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in  theory.  However,  it  turned  out  to  be  very  fast  and  accurate  and  worked  in  all  of  the  practical 
range  of  physically-relevant  cases  (besides,  it  can  easily  be  modified  to  be  made  systematic).  It 
consists  in  choosing  an  initial  step  size  and  computing  the  error  function  along  consecutive 

values  within  the  segment [^0  -d<pld£/,£0  +  d(/)ldCr\,  until  we  find  three  consecutive  values 

that  C-i  ’Cv’Cv+i  such  that  <b(C-i)>(^(C)  and  ^>(Cv)<^>(C+i)-  We  theh  reduce  the  step  size 
and  continue  until  the  segment  bracketing  the  minimum  is  as  narrow  as  desired. 

Gradient  computation  and  line  minimization  are  then  performed  repeatedly  until  the  gradient 
becomes  null  (and  we  have  reached  the  optimal  source  distribution)  or,  most  likely,  until  the  last 
iterations  results  in  an  insignificant  change  of  the  objective  function  (i.e.  below  a  given 
threshold) 

Key  Research  Accomplishment : 

For  a  total  target  size  of  50*50,  which  is  already  quite  simplistic,  the  typical  computation  time 
for  about  10  iterations  was  initially  of  6  hours  (on  Matlab,  with  a  400MHz  PC)  and  then 
optimized  down  to  2  hours,  which  is  still  cumbersome. 

We  started  with  a  single  light-emitting  source  located  in  a  2D  array  at  (0,0)  (i.e.  in  the  center  of 
50*50  matrix).  The  reconstruction  scheme  was  initialized  from  a  distribution  of  9  sources  located 
in  the  square  [-1;1]*[-1;1].  For  each  step,  we  assumed  that  we  could  only  access  the  intensity  of 
light  emitted  at  the  boundary  of  the  matrix,  just  like  a  photon-counting  device  would  only  be  able 
to  count  the  photons  at  the  boundary  of  the  medium 

The  source  distribution  obtained  after  100  iterations  by  the  method  is  shown  in  figure  2. 
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Source  Distribution  Reconstructed  for  ja  Single. 
Source  in  (0,0);  starting  from  a  distribution,  j 
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fig  2. :  This  is  the  source  distribution  that  was  reconstructed  for 
a  single  source  located  at  the  center  of  the  matrix  (the  x  and  y  axes  only  reflect  the  indices  of  the 
matrix  elements,  not  their  ‘actual’  location),  departing  from  an  initial  source  distribution  that  was 
located  in  the  3*3  voxels  around  the  center  of  the  matrix. 


Conclusion: 

The  basic  iteration  scheme  can  be  applied  to  this  problem.  Further  work  is  needed  to  improve  the 
forward  calculation  to  make  it  more  physically  realistic,  and,  in  particular,  the  finite-difference 
scheme  needs  to  be  incorporated  in  the  gradient  computation.  Future  work  should  include  the 
investigation  of  algorithm  when  the  optical  attenuation  jua ,  optical  scatter  jus  are  spatially 
dependent. 
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